

% who are you?
% analyst = 1,...,4;


%%
load(['sim_eq_round1_',num2str(analyst),'a_replication.mat'])
rng(srng)

% strategy profiles
YY=zeros(1,6);
for i5=0:1
    for i4=0:1
        for i3=0:1
            for i2=0:1
                for i1=0:1
                    YY=[YY;[1,i1,i2,i3,i4,i5]]; %#ok<AGROW>
                end
            end
        end
    end
end
clear i1 i2 i3 i4 i5


%% set mean utilities
load('data.mat')
X=ones(size(data,1),1); [N,K]=size(X);
Z_US=zeros(N,1); KZ=size(Z_US,2);
Z_UK=zeros(N,1);
Z_France=zeros(N,1);
Z_Russia=zeros(N,1);
Z_China=zeros(N,1);
Z=[Z_US,Z_UK,Z_France,Z_Russia,Z_China];
reb={1:K,[]}; mpow={1,1:KZ};
thetaStar = zeros(length(reb{1})+4+length(mpow{1})+5+length(reb{2})+length(mpow{2})+10,1);
U0 = payoffs(thetaStar,X,Z,YY,reb,mpow);
clearvars -except analyst N srng U0 YY


%% simulated equilibria 
N1 = floor(N/4);
NS = 1e+3;

if analyst<4
    
    V_R = U0(:,1,N1*(analyst-1)+(1:N1)) + randn(size(YY,1),NS,N1); % entry payoffs of rebel group
    V_US = U0(:,2,N1*(analyst-1)+(1:N1)) + randn(size(YY,1),NS,N1);
    V_UK = U0(:,3,N1*(analyst-1)+(1:N1)) + randn(size(YY,1),NS,N1);
    V_France = U0(:,4,N1*(analyst-1)+(1:N1)) + randn(size(YY,1),NS,N1);
    V_Russia = U0(:,5,N1*(analyst-1)+(1:N1)) + randn(size(YY,1),NS,N1);
    V_China = U0(:,6,N1*(analyst-1)+(1:N1)) + randn(size(YY,1),NS,N1);

    filen=1:NS;
    intervEq=cell(NS,N1);

    for n=1:N1
        tic
        parfor b=1:NS
            fn=filen(b);
            fname=['game',num2str(analyst),'a',num2str(fn),'.nfg'];
            fileID=fopen(fname,'w');
            fprintf(fileID,'NFG 1 R "Major-Powers Intervention"\n',1); %#ok<CTPCT>
            fprintf(fileID,'{ "US" "UK" "France" "Russia" "China" } { 2 2 2 2 2 }\n',1); %#ok<CTPCT>
            formatSpec='%12g ';
            fprintf(fileID,formatSpec,reshape([V_US(2:end,b,n),V_UK(2:end,b,n),V_France(2:end,b,n),V_Russia(2:end,b,n),V_China(2:end,b,n)]',1,5*2^5)); %#ok<*PFBNS>
            fclose(fileID);        
            [~,NE]=system(['gambit-enumpoly ',fname]);
            findNE=strfind(NE,'NE');
            EQ=zeros(2^5,length(findNE));
            EV_R=zeros(1,length(findNE));
            for m=1:length(findNE)
                sigma_US=str2double(NE(findNE(m)+12:findNE(m)+19));
                sigma_UK=str2double(NE(findNE(m)+30:findNE(m)+37));
                sigma_France=str2double(NE(findNE(m)+48:findNE(m)+55));
                sigma_Russia=str2double(NE(findNE(m)+66:findNE(m)+73));
                sigma_China=str2double(NE(findNE(m)+84:findNE(m)+91));
                EQ(:,m)=(sigma_US.^YY(2:end,2)).*((1-sigma_US).^(1-YY(2:end,2))).*...
                    (sigma_UK.^YY(2:end,3)).*((1-sigma_UK).^(1-YY(2:end,3))).*...
                    (sigma_France.^YY(2:end,4)).*((1-sigma_France).^(1-YY(2:end,4))).*...
                    (sigma_Russia.^YY(2:end,5)).*((1-sigma_Russia).^(1-YY(2:end,5))).*...
                    (sigma_China.^YY(2:end,6)).*((1-sigma_China).^(1-YY(2:end,6)));
                EV_R(m)=V_R(2:end,b,n)'*EQ(:,m);
            end
            intervEq{b,n}=[EV_R;EQ];
        end
        toc
    end

else
    
    V_R = U0(:,1,3*N1+1:N) + randn(size(YY,1),NS,N-3*N1); % entry payoffs of rebel group
    V_US = U0(:,2,3*N1+1:N) + randn(size(YY,1),NS,N-3*N1);
    V_UK = U0(:,3,3*N1+1:N) + randn(size(YY,1),NS,N-3*N1);
    V_France = U0(:,4,3*N1+1:N) + randn(size(YY,1),NS,N-3*N1);
    V_Russia = U0(:,5,3*N1+1:N) + randn(size(YY,1),NS,N-3*N1);
    V_China = U0(:,6,3*N1+1:N) + randn(size(YY,1),NS,N-3*N1);

    filen=1:NS;
    intervEq=cell(NS,N-3*N1);

    for n=1:N-3*N1
        tic
        parfor b=1:NS
            fn=filen(b);
            fname=['game',num2str(analyst),'a',num2str(fn),'.nfg'];
            fileID=fopen(fname,'w');
            fprintf(fileID,'NFG 1 R "Major-Powers Intervention"\n',1); %#ok<CTPCT>
            fprintf(fileID,'{ "US" "UK" "France" "Russia" "China" } { 2 2 2 2 2 }\n',1); %#ok<CTPCT>
            formatSpec='%12g ';
            fprintf(fileID,formatSpec,reshape([V_US(2:end,b,n),V_UK(2:end,b,n),V_France(2:end,b,n),V_Russia(2:end,b,n),V_China(2:end,b,n)]',1,5*2^5));
            fclose(fileID);        
            [~,NE]=system(['gambit-enumpoly ',fname]);
            findNE=strfind(NE,'NE');
            EQ=zeros(2^5,length(findNE));
            EV_R=zeros(1,length(findNE));
            for m=1:length(findNE)
                sigma_US=str2double(NE(findNE(m)+12:findNE(m)+19));
                sigma_UK=str2double(NE(findNE(m)+30:findNE(m)+37));
                sigma_France=str2double(NE(findNE(m)+48:findNE(m)+55));
                sigma_Russia=str2double(NE(findNE(m)+66:findNE(m)+73));
                sigma_China=str2double(NE(findNE(m)+84:findNE(m)+91));
                EQ(:,m)=(sigma_US.^YY(2:end,2)).*((1-sigma_US).^(1-YY(2:end,2))).*...
                    (sigma_UK.^YY(2:end,3)).*((1-sigma_UK).^(1-YY(2:end,3))).*...
                    (sigma_France.^YY(2:end,4)).*((1-sigma_France).^(1-YY(2:end,4))).*...
                    (sigma_Russia.^YY(2:end,5)).*((1-sigma_Russia).^(1-YY(2:end,5))).*...
                    (sigma_China.^YY(2:end,6)).*((1-sigma_China).^(1-YY(2:end,6)));
                EV_R(m)=V_R(2:end,b,n)'*EQ(:,m);
            end
            intervEq{b,n}=[EV_R;EQ];
        end
        toc
    end

end


%%
clear filen n
save(['sim_eq_round1_',num2str(analyst),'a.mat'])


